
FUNCTION chisq_tnmin_f, A0, XARR=xarr, YARR=yarr, SIGHIST=sighist, N=n, $
                        MU=mu, SIGMA=sig, _EXTRA=_extra
  
  ln = gauss_1(ALOG(10.^xarr),[A0,mu,sig])
  chisq = total( (yarr-ln)*(yarr-ln)/sighist/sighist )
  print,chisq / (n-2)
  
  RETURN, chisq
END


PRO chisq_ln_fit
  
  mw       = omni_read_conffile('./conffiles/galactic_params.conf')
  
  mu  = [4.15,5.20]
  sig = [1.36,1.19]
  tdstr = [string(mw.td,format="(I0)")+' K','Log-Normal']
  mmin = [108.5,50.4]
  
  
  FOR td=0,1 DO BEGIN
     
     
     restore,'./masses/mfn_mega_td'+strtrim(td,2)+'.sav',/ver
     
     
     
     ;; Do histogramming here...
     nreal  = n_elements(dist)
     ndloop = nreal
     mconf  = omni_read_conffile('./conffiles/mass_deriv.conf')
     
     ;; Plot the Differential Mass Function
     binw = mconf.dmfbinw NE 0. ? mconf.dmfbinw : fdr_width(alog10(mfn))
     plothist,alog10(mfn),xarr,yarr,bin=binw,/noplot
     yarr = [0,yarr,0] / float(ndloop)
     dx   = mean(xarr[1:*]-xarr[0,*])
     xarr = [xarr[0]-dx,xarr,xarr[n_elements(xarr)-1]+dx]
     
     mmin[td] = mmin[td] > 10.^(min(xarr)+dx)
     
     message,'Computing differential mass function error bars...',/inf
     nbin = n_elements(xarr)
     hist = fltarr(nreal,nbin)
     ;; For each realization, make distance cuts & histogram
     FOR ii=0, nreal-1 DO BEGIN
        d = dist[ii].objs / 1.d3
        dind = where(d GE mconf.dmin AND d LE mconf.dmax, nd)
        hist[ii,*] = histogram(alog10(dist[ii].mass[dind]), binsize=binw,$
                               min=min(xarr)-dx/2., nbins=nbin, loc=xh)
     ENDFOR
     ;; Loop through bins to get mu, sig
     sighist = fltarr(nbin)
     FOR jj=0,nbin-1 DO sighist[jj] = stddev(hist[*,jj])
     
     ;; Go go plotting fun!
     myps,'./masses/plots/ln_chisq_td'+strtrim(td,2)+'.eps'
     
     sun = cgSymbol('sun')
     cgPlot,ALOG(10.^xarr),yarr,psym=10,/ylog,yr=[1d-1,1d2],/yst,$
            ymargin=[4,4],ytickformat='exponent10',xtit='ln M'+sun,thick=3,$
            charsize=1.0,xst=9,ytit='Mean N per log(M) = 0.1 bin',xr=[0,10.5]
     oploterror,ALOG(10.^xarr),yarr,sighist,psym=3,errthick=2
     
     oploterror,ALOG(10.^xarr),yarr,sqrt(yarr),psym=3,color='cyan',errthick=2
     
     cgAxis,xaxis=1,xr=exp(!x.crange),/xst,xtit='M  [M'+sun+']',/xlog,$
            charsize=1.0,xtickformat='exponent10'
     
     vline,alog(mmin[td]),/log
     vline,alog(8d3),/log
     
     
     ind = where(10.^xarr GE mmin[td] AND 10.^xarr LE 8000, n)
     print,n
     
     yarr    = yarr[ind]
     xarr    = xarr[ind]
     sighist = sighist[ind]
     
     A0 = 30.
     
     fcnargs = {XARR:xarr,$
                YARR:yarr,$
                SIGHIST:sighist,$
                N:n,$
                MU:mu[td],$
                SIGMA:sig[td]}
     
     norm = TNMIN('chisq_tnmin_f',A0,FUNCTARGS=fcnargs,STATUS=status,$
                  ERRMSG=errmsg,/AUTODERIVATIVE,BESTMIN=bestmin)
     
     print,norm
     print,status
     print,errmsg
     
     ln = gauss_1(ALOG(10.^xarr),[norm,mu[td],sig[td]])
     cgOplot,ALOG(10.^xarr),ln,color='red',thick=5
     
     sighist_sdm = sighist / sqrt(nreal)
     chisq_sdm = total( (yarr-ln)*(yarr-ln)/ sighist_sdm / sighist_sdm )
     
     print,'CHISQ_SDM: ',chisq_sdm / (n-2)
     
     al_legend,/top,/left,/clear,$
               [cgSymbol('chi')+'!u2!dred!n = '+$
                string(bestmin/(n-2),format="(F0.2)"),$
                'T!dd!n = '+tdstr[td],$
                cgSymbol('mu')+' = '+string(mu[td],format="(F0.2)"),$
                cgSymbol('sigma')+' = '+string(sig[td],format="(F0.2)"),$
                'N!dreal!n = '+string(nreal,format="(I0)")]
     
     al_legend,/top,/right,color=['black','cyan'],thick=2,linestyle=0,$
               ['Dispersion','Poisson'],linsize=0.5,charsize=0.8,/clear
     
     myps,/done
     
  ENDFOR
  
END
  
